This invention relates to radar systems, and especially to schemes for reducing the noise in echoes returned from diffuse targets such as weather phenomena.
The use of radar to detect major weather phenomena such as thunderstorms is well known to the public at large through television weather reports. Weather researchers and air traffic controllers use radar to detect and track low-amplitude or low level echoes from disturbances less severe than thunderstorms or precipitation. Such low level disturbances include small changes in the dielectric constant or index of refraction of disturbed atmosphere, believed to be attributable to slight temperature, humidity, moisture or pressure variations, but which may also be attributable to the presence of insects, or dust associated with the disturbance.
Copending patent application Ser. No. 07/685,792, filed Apr. 16, 1991 in the name of Urkowitz, describes an air traffic control radar which is also adapted for detecting low-amplitude weather phenomena. As described therein, a pulse doppler radar includes range sidelobe reduction following the doppler filtering, in order to reduce range errors. Another arrangement for reducing range sidelobes in a pulse doppler radar, including the transmission of complementary-sequence pulses, with matched filtering following the doppler filtering, is described in copending application Ser. No. 07/734,003 filed Jul. 22, 1991 in the name of Urkowitz.
FIG. 1 is a simplified block diagram of a portion of a radar system as generally described in the abovementioned Urkowitz patent applications. In its simplified form, it may be taken as corresponding to the prior art. In general, radar system 10 of FIG. 1 includes a transmit-receive multiplexer (TR) apparatus 14. A transmitter (TX) 16 and a receiver 18 are connected to antenna 12 by way of TR 14. Receiver 18 is connected to transmitter 16 for receiving reference local oscillator (LO) signals. An analog-to-digital converter (ADC) 20 is connected to the output of receiver 18, and is in turn connected to a digital signal processor 22. In operation, transmitter 16 produces sequences of pulses which are transmitted by antenna 12 in the form of electromagnetic radiation or signal, represented as 24. The radiation reaches a weather phenomenon designated generally by 26, which reflects signal echoes back toward radar system 10 and its antenna 12. Phenomenon 26 may be turbulent, with internal fluid flow in whorls, illustrated by arrows 25. Antenna 12 couples the received echo signals by way of TR 14 to receiver 18. Receiver 18 also receives local oscillator (LO) signal from transmitter 16 to aid in its processing. The received signals are sampled and converted into digital form (parallel or serial) in ADC 20. The digital signals are applied to a digital signal processor (DSP) 22 for processing the signals to extract information therefrom. Many of the processes which are performed in DSP 22 are performed by algorithms rather than by dedicated hardware in modern radar systems, but the operation of the algorithms is often described in terms of their analog hardware equivalents. Thus, for example, the discrete Fourier transform (DFT) is often used to identify the frequency components of a series of pulses. This process may be described by its analog equivalent, which is applying the analog signal to a bank of narrow-band filters. The signals processed in DSP 22 by DFT and by other processes are then made available for further processing, such as generating a surveillance display, activating alarms, and the like. Such additional processing steps are well known in the art and are not elaborated herein.
The extremely low amplitudes of the echoes from weather phenomena creates difficulties in identifying returns from the targets and in distinguishing them from atmospheric and radar system noise. If the amplitude of the atmospheric and system noise were known beforehand, a fixed amount could be subtracted from each signal sample, to thereby effectively render the radar noise-free. However, the noise level of the radar system (the transmitter and the receiver) itself changes from time to time, and atmospheric noise effects may change drastically from moment to moment. Such changes may occur, for example, as a result of arc discharges such as automobile spark plugs or electric motors, which produce broadband impulse noise.
An article entitled "Objective Determination of the Noise Level in Doppler Spectra", by Hildebrand and Sekhon, which was published in the Journal of Applied Meteorology, October 1974, pages 808-811 (the HS article or HS), describes a noise reduction scheme which may be represented by FIG. 2. FIG. 2 represents one of the first stages of processing in DSP 22 of the radar of FIG. 1. In FIG. 2, digitized echoes resulting from sequences of transmitted pulses are applied to a doppler filter bank 210, including filters f.sub.1, f.sub.2. . . f.sub.n, which separate signals representative of the echoes or returns from a weather phenomenon or other target into their various frequency components f.sub.1, f.sub.2, . . . f.sub.n. These doppler frequency components represent echoes from portions of the diffuse target which are moving with different radial velocities relative to the radar system. If a single solid object such as a cannon-ball were moving with a fixed radial velocity relative to the radar, only one output signal, representing one radial velocity, would be expected from doppler filter bank 210. However, a diffuse target with turbulence may have different portions moving with different velocities, and these portions contribute to the overall doppler signal at the output of filter 210.
In FIG. 2, the output of each filter f.sub.1, f.sub.2, . . . f.sub.n of doppler filter bank 210 is applied to the noninverting (+) input of a corresponding summing circuit 212. Thus, filter portion f.sub.1 of doppler filter bank 210 is coupled to the noninverting input of a summing circuit 212.sup.1, filter f.sub.2 is coupled to a summing circuit 212.sup.2, . . . and the n.sup.th filter f.sub.n is coupled to a summing circuit 212.sub.n. A noise determination or noise calculation circuit 214 is connected to the output of each doppler filter of filter bank 210 and to the inverting (-) input of the corresponding summing circuit 212. Thus, a noise calculation circuit 214 is connected to the output of filters f.sub.1, f.sub.2, . . . f.sub.n, and calculates the overall noise in channel 1, and applies the calculated noise in common to the inverting inputs of the summing circuits 212.sup.1. . . 212.sup.n. Within each summing circuit 212, the calculated noise is subtracted from the signal-plus-noise applied to its noninverting input, to thereby produce ideally noise-free signal.
The Hildebrand, et al. (HS) article describes a method of determining the noise level which depends upon the physical properties of white and Gaussian noise. The HS method for determining the noise level is described in conjunction with FIGS. 3a-3h.
The output signal from ADC 20 of FIG. 1 is a sequence of pulse-to-pulse complex values from each range bin. FIG. 3a plots received signal amplitude x(f) versus frequency f at the outputs of the doppler filters. Plot 310 represents a possible signal amplitude spectrum within a single range bin resulting from transmission of a sequence of pulses. Curve 310 of FIG. 3a, when squared, represents the power density spectrum S(f) of the received signal. FIG. 3b plots S(f), the power density spectrum of the signal across the doppler channels of FIG. 2. The value of S(f) at any point may be given by EQU S(f).varies.e.sup.-(f-f.sub.d.sup.).sup.2.sub./2.sigma..sup.2.sub.f ( 1)
where .varies. represents proportionality, f is the frequency f.sub.d is the doppler frequency and .sigma..sub.f is the spread or standard deviation. Curve 312 is generally bell-shaped. Within any frequency increment f.sub.i of FIG. 3b, the mean-square value of the power density spectrum represents the power within that frequency increment, and the sum of the incremental powers represented by all such frequency increments represents the signal power in the doppler channels. The spread or standard deviation of curve 312 of FIG. 3b is represented by arrow .sigma..sub.f. Ordinarily, each frequency increment f.sub.i corresponds to the bandwidth of one doppler filter.
The signal component of FIG. 3b is accompanied by a constant noise power density spectrum N.sub.o, illustrated in FIG. 3c by plot 314. The magnitude of N.sub.o is not known, but is expected to vary with time. Plot 314 of FIG. 3c represents the constant noise power density spectrum across all the doppler filter channels having a cumulative bandwidth F. The spectral variance of the noise spectrum of FIG. 3b is given by EQU .sigma..sub.N.sup.2 =F.sup.2 /12 (2)
where F is the pulse recurrence frequency (PRF). Equation 2 corresponds to the moment of inertia of the rectangle defined by plot 314 about center frequency F/2 of FIG. 3b. Plot 316 of FIG. 3d represents the sum of the power density spectrums of FIGS. 3b and 3c, which is the echo signal plus noise within one range bin resulting from a sequence of pulses.
As a result of the use of digital signal processing, the signals represented in FIGS. 3a, 3b, 3c and 3d are actually time and amplitude-quantized. Further, the discrete Fourier transform (DFT) by which the doppler filtering is actually accomplished produces discrete sample points which generally number 2.sup.N. For example, the number of DFT sample points might be 16, 32. . . . FIG. 3e plots a plurality of points which might be produced by a DFT in response to signal equivalent to 316 of FIG. 3d. For simplicity, only 8 points, numbered from 1 to 8, are represented in FIG. 3d, rather than the larger number which would ordinarily be used. These points represent samples of analog plot 316 of FIG. 3d, taken at various frequencies.
In the prior-art HS method for calculating noise in each channel, the samples represented as samples 1-8 in FIG. 3e are ordered according to amplitude, or in effect laid out in a line, with smallest-amplitude samples 1 and 8 at one end, and largest-amplitude sample 5 at the other end. FIG. 3f represents such an ordering of samples, with highest amplitudes near the top and lowest amplitude near the bottom. For the purpose of calculating the noise level, one or more of the largest-amplitude samples is discarded. Since only eight samples are illustrated in the simplified example, only one sample is deleted. The largest-amplitude sample, which is sample 5, is deleted to make a reduced-order set of samples (periodogram). This is the equivalent of procrustinating the quantized plot of FIG. 3e by truncating the plot at a level illustrated by dashed line 318, to produce a plot which, if not quantized, might look like plot 322 of FIG. 3g.
When one or more samples of the quantized plot of FIG. 3e, such as sample 5, are discarded, the remaining samples have a "frequency gap", which happens to be near fd in FIG. 3e. According to the HS method for calculating the noise, the frequency gap is "closed" to form a continuous spectrum, illustrated as 324 in FIG. 3h. A test is then performed to see if the remaining spectrum has the characteristics attributed to white noise. This test involves the calculation of two ratios, R1 and R2, where EQU R.sub.1 =.sigma..sup.2.sub.N /.sigma..sup.2 ( 3) EQU R.sub.2 =P.sup.2 /pQ (4)
The quantity .sigma..sub.N.sup.2, the presumed noise spectral variance, is given by equation (2). The other quantities are expressed as EQU .sigma..sup.2 =(.SIGMA.f.sup.2.sub.n S.sub.n /.SIGMA.S.sub.n)-(.SIGMA.f.sub.n S.sub.n /.SIGMA.S.sub.n).sup.2 ( 5) EQU P=.SIGMA.S.sub.n /N (6) EQU Q=.SIGMA.(S.sup.2.sub.n /N)-P.sup.2 ( 7)
where
f.sub.n is the frequency of the specific one of the N doppler filter outputs PA1 S.sub.n is the value of estimated spectral density at frequency f.sub.n PA1 F is the frequency spread of the spectrum, PA1 N is the number of independent spectral densities or doppler filters in the filter bank, and PA1 p is the number of spectral lines over which a moving average is taken.
Ratios R.sub.1 and R.sub.2 are calculated after the reduction in the order of the set of samples. If the HS test is not met, one or more of the largest amplitude samples are deleted from the reduced-order set, to produce a new reduced-order set, and the test is applied again. The procedure of reducing the order and performing the test is performed as many times as necessary. According to the HS test, the remaining samples represent white noise when both ratios R.sub.1 and R.sub.2 equal unity.
When R.sub.1 and R.sub.2 both sufficiently close to unity, the samples remaining in the current reduced-order sample set are taken to represent the noise level. The resulting noise spectral density level is subtracted from the received samples in each channel, to produce what is expected to be a more accurate representation of the echo signal. This calculation and subtraction is performed for each range bin.
The HS scheme suffers from both theoretical and practical problems. The theoretical problem is that the HS calculation assumes that the noise is Gaussian. This assumption may not be true, especially when the noise originates from non-thermal phenomenon. Non-thermal phenomena might include phase noise, possibly attributable to timing jitter, or it might include impulse noise such as automobile ignition noise, lightning, electric motors or the like. Also, the ratios R.sub.1 and R.sub.2 are unlikely to reach exactly unity at the same time. The practical problem which then arises is how to determine how close each must be to unity in order for the condition to be deemed to be fulfilled, and the determination of the confidence level for whatever deviation from unity is selected. Thus, the HS test provides no standards by which a time for completion of the iterations can be determined.